{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "using VMLS\n",
    "using LinearAlgebra"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Chapter 9\n",
    "# Linear dynamical systems\n",
    "### 9.1 Linear dynamical systems\n",
    "Let’s simulate a time-invariant linear dynamic system \n",
    "$$\n",
    "x_{t+1} = Ax_t, t = 1, . . . , T,\n",
    "$$\n",
    "with dynamics matrix\n",
    "$$\n",
    "A =\n",
    "\\begin{bmatrix}\n",
    "0.97 &  0.10  & −0.05 \\\\\n",
    "−0.3 &  0.99  & 0.05 \\\\\n",
    "0.01 &  −0.04 &  0.96\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "and initial state $x_1 = (1, 0,−1)$. We store the state trajectory in the $n × T$ matrix\n",
    "`state_traj`, with the $i$th column $x_t$. We plot the result in Figure 9.1."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/svg+xml": [
       "<?xml version=\"1.0\" encoding=\"utf-8\"?>\n",
       "<svg xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" width=\"600\" height=\"400\" viewBox=\"0 0 2400 1600\">\n",
       "<defs>\n",
       "  <clipPath id=\"clip5000\">\n",
       "    <rect x=\"0\" y=\"0\" width=\"2400\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip5000)\" d=\"\n",
       "M0 1600 L2400 1600 L2400 0 L0 0  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip5001\">\n",
       "    <rect x=\"480\" y=\"0\" width=\"1681\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip5000)\" d=\"\n",
       "M141.865 1425.62 L2352.76 1425.62 L2352.76 47.2441 L141.865 47.2441  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip5002\">\n",
       "    <rect x=\"141\" y=\"47\" width=\"2212\" height=\"1379\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  161.871,1425.62 161.871,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  587.533,1425.62 587.533,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1013.2,1425.62 1013.2,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1438.86,1425.62 1438.86,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1864.52,1425.62 1864.52,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  2290.18,1425.62 2290.18,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  141.865,1266.26 2352.76,1266.26 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  141.865,984.236 2352.76,984.236 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  141.865,702.213 2352.76,702.213 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  141.865,420.19 2352.76,420.19 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  141.865,138.167 2352.76,138.167 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  141.865,1425.62 2352.76,1425.62 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  141.865,1425.62 141.865,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  161.871,1425.62 161.871,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  587.533,1425.62 587.533,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1013.2,1425.62 1013.2,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1438.86,1425.62 1438.86,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1864.52,1425.62 1864.52,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  2290.18,1425.62 2290.18,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  141.865,1266.26 175.028,1266.26 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  141.865,984.236 175.028,984.236 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  141.865,702.213 175.028,702.213 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  141.865,420.19 175.028,420.19 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  141.865,138.167 175.028,138.167 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 161.871, 1479.62)\" x=\"161.871\" y=\"1479.62\">0</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 587.533, 1479.62)\" x=\"587.533\" y=\"1479.62\">10</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1013.2, 1479.62)\" x=\"1013.2\" y=\"1479.62\">20</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1438.86, 1479.62)\" x=\"1438.86\" y=\"1479.62\">30</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1864.52, 1479.62)\" x=\"1864.52\" y=\"1479.62\">40</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 2290.18, 1479.62)\" x=\"2290.18\" y=\"1479.62\">50</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 117.865, 1283.76)\" x=\"117.865\" y=\"1283.76\">-2</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 117.865, 1001.74)\" x=\"117.865\" y=\"1001.74\">-1</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 117.865, 719.713)\" x=\"117.865\" y=\"719.713\">0</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 117.865, 437.69)\" x=\"117.865\" y=\"437.69\">1</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 117.865, 155.667)\" x=\"117.865\" y=\"155.667\">2</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:66px; text-anchor:middle;\" transform=\"rotate(0, 1247.31, 1559.48)\" x=\"1247.31\" y=\"1559.48\">t</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#009af9; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  204.437,420.19 247.003,414.55 289.569,419.654 332.136,435.354 374.702,461.15 417.268,496.212 459.834,539.399 502.401,589.294 544.967,644.252 587.533,702.449 \n",
       "  630.1,761.942 672.666,820.73 715.232,876.823 757.798,928.306 800.365,973.401 842.931,1010.53 885.497,1038.36 928.063,1055.87 970.63,1062.37 1013.2,1057.51 \n",
       "  1055.76,1041.33 1098.33,1014.26 1140.89,977.06 1183.46,930.866 1226.03,877.101 1268.59,817.46 1311.16,753.842 1353.73,688.3 1396.29,622.964 1438.86,559.979 \n",
       "  1481.42,501.429 1523.99,449.271 1566.56,405.268 1609.12,370.929 1651.69,347.455 1694.26,335.702 1736.82,336.145 1779.39,348.86 1821.95,373.52 1864.52,409.402 \n",
       "  1907.09,455.403 1949.65,510.082 1992.22,571.696 2034.79,638.261 2077.35,707.612 2119.92,777.474 2162.48,845.539 2205.05,909.536 2247.62,967.312 2290.18,1016.9 \n",
       "  \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#e26f46; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  204.437,702.213 247.003,800.921 289.569,899.629 332.136,994.942 374.702,1083.55 417.268,1162.36 459.834,1228.57 502.401,1279.78 544.967,1314.06 587.533,1330.06 \n",
       "  630.1,1326.99 672.666,1304.73 715.232,1263.75 757.798,1205.18 800.365,1130.75 842.931,1042.69 885.497,943.759 928.063,837.058 970.63,725.994 1013.2,614.144 \n",
       "  1055.76,505.146 1098.33,402.576 1140.89,309.83 1183.46,230.01 1226.03,165.821 1268.59,119.475 1311.16,92.6139 1353.73,86.2547 1396.29,100.745 1438.86,135.751 \n",
       "  1481.42,190.256 1523.99,262.597 1566.56,350.505 1609.12,451.183 1651.69,561.387 1694.26,677.536 1736.82,795.821 1779.39,912.333 1821.95,1023.19 1864.52,1124.67 \n",
       "  1907.09,1213.31 1949.65,1286.08 1992.22,1340.42 2034.79,1374.36 2077.35,1386.61 2119.92,1376.56 2162.48,1344.35 2205.05,1290.85 2247.62,1217.63 2290.18,1126.95 \n",
       "  \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5002)\" style=\"stroke:#3da44d; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  204.437,984.236 247.003,970.135 289.569,952.593 332.136,931.856 374.702,908.292 417.268,882.385 459.834,854.712 502.401,825.93 544.967,796.749 587.533,767.914 \n",
       "  630.1,740.175 672.666,714.262 715.232,690.865 757.798,670.603 800.365,654.01 842.931,641.509 885.497,633.401 928.063,629.853 970.63,630.89 1013.2,636.393 \n",
       "  1055.76,646.102 1098.33,659.62 1140.89,676.43 1183.46,695.905 1226.03,717.332 1268.59,739.932 1311.16,762.885 1353.73,785.358 1396.29,806.532 1438.86,825.625 \n",
       "  1481.42,841.925 1523.99,854.807 1566.56,863.758 1609.12,868.396 1651.69,868.477 1694.26,863.912 1736.82,854.766 1779.39,841.258 1821.95,823.758 1864.52,802.77 \n",
       "  1907.09,778.922 1949.65,752.942 1992.22,725.636 2034.79,697.866 2077.35,670.515 2119.92,644.461 2162.48,620.55 2205.05,599.564 2247.62,582.198 2290.18,569.033 \n",
       "  \n",
       "  \"/>\n",
       "<path clip-path=\"url(#clip5000)\" d=\"\n",
       "M1890.96 372.684 L2280.76 372.684 L2280.76 130.764 L1890.96 130.764  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1890.96,372.684 2280.76,372.684 2280.76,130.764 1890.96,130.764 1890.96,372.684 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#009af9; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1914.96,191.244 2058.96,191.244 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:start;\" transform=\"rotate(0, 2082.96, 208.744)\" x=\"2082.96\" y=\"208.744\">(x_t)_1</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#e26f46; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1914.96,251.724 2058.96,251.724 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:start;\" transform=\"rotate(0, 2082.96, 269.224)\" x=\"2082.96\" y=\"269.224\">(x_t)_2</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip5000)\" style=\"stroke:#3da44d; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1914.96,312.204 2058.96,312.204 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:start;\" transform=\"rotate(0, 2082.96, 329.704)\" x=\"2082.96\" y=\"329.704\">(x_t)_3</text>\n",
       "</g>\n",
       "</svg>\n"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x_1 = [1,0,-1]; # initial state\n",
    "n = length(x_1); T = 50;\n",
    "A = [ 0.97 0.10 -0.05 ; -0.3 0.99 0.05 ; 0.01 -0.04 0.96 ]\n",
    "state_traj = [x_1 zeros(n,T-1) ];\n",
    "for t=1:T-1 # Dynamics recursion\n",
    "state_traj[:,t+1] = A*state_traj[:,t];\n",
    "end\n",
    "using Plots\n",
    "plot(1:T, state_traj', xlabel = \"t\", label = [\"(x_t)_1\", \"(x_t)_2\", \"(x_t)_3\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Figure 9.1** Linear dynamical system simulation.\n",
    "\n",
    "### 9.2 Population dynamics\n",
    "We can create a population dynamics matrix with just one simple line of Julia. The following code predicts the 2020 population distribution in the US using the data of Section [9.2](https://web.stanford.edu/~boyd/vmls/vmls.pdf#section.9.2) of VMLS, which are available through the VMLS function `population_data`. The result is shown in Figure 9.2."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/svg+xml": [
       "<?xml version=\"1.0\" encoding=\"utf-8\"?>\n",
       "<svg xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" width=\"600\" height=\"400\" viewBox=\"0 0 2400 1600\">\n",
       "<defs>\n",
       "  <clipPath id=\"clip5400\">\n",
       "    <rect x=\"0\" y=\"0\" width=\"2400\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip5400)\" d=\"\n",
       "M0 1600 L2400 1600 L2400 0 L0 0  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip5401\">\n",
       "    <rect x=\"480\" y=\"0\" width=\"1681\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip5400)\" d=\"\n",
       "M175.611 1425.62 L2352.76 1425.62 L2352.76 47.2441 L175.611 47.2441  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip5402\">\n",
       "    <rect x=\"175\" y=\"47\" width=\"2178\" height=\"1379\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  216.481,1425.62 216.481,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  735.146,1425.62 735.146,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1253.81,1425.62 1253.81,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1772.47,1425.62 1772.47,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  2291.14,1425.62 2291.14,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  175.611,1398.89 2352.76,1398.89 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  175.611,1110.08 2352.76,1110.08 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  175.611,821.268 2352.76,821.268 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  175.611,532.455 2352.76,532.455 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  175.611,243.643 2352.76,243.643 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  175.611,1425.62 2352.76,1425.62 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  175.611,1425.62 175.611,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  216.481,1425.62 216.481,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  735.146,1425.62 735.146,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1253.81,1425.62 1253.81,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1772.47,1425.62 1772.47,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  2291.14,1425.62 2291.14,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  175.611,1398.89 208.268,1398.89 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  175.611,1110.08 208.268,1110.08 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  175.611,821.268 208.268,821.268 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  175.611,532.455 208.268,532.455 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5400)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  175.611,243.643 208.268,243.643 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 216.481, 1479.62)\" x=\"216.481\" y=\"1479.62\">0</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 735.146, 1479.62)\" x=\"735.146\" y=\"1479.62\">25</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1253.81, 1479.62)\" x=\"1253.81\" y=\"1479.62\">50</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1772.47, 1479.62)\" x=\"1772.47\" y=\"1479.62\">75</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 2291.14, 1479.62)\" x=\"2291.14\" y=\"1479.62\">100</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 151.611, 1416.39)\" x=\"151.611\" y=\"1416.39\">0</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 151.611, 1127.58)\" x=\"151.611\" y=\"1127.58\">1</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 151.611, 838.768)\" x=\"151.611\" y=\"838.768\">2</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 151.611, 549.955)\" x=\"151.611\" y=\"549.955\">3</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 151.611, 261.143)\" x=\"151.611\" y=\"261.143\">4</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:66px; text-anchor:middle;\" transform=\"rotate(0, 1264.18, 1559.48)\" x=\"1264.18\" y=\"1559.48\">Age</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5400)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:66px; text-anchor:middle;\" transform=\"rotate(-90, 89.2861, 736.431)\" x=\"89.2861\" y=\"736.431\">Population (millions)</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip5402)\" style=\"stroke:#009af9; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  237.228,220.804 257.974,228.178 278.721,227.836 299.468,228.322 320.214,229.489 340.961,229.997 361.707,233.913 382.454,237.236 403.2,241.111 423.947,244.925 \n",
       "  444.694,268.725 465.44,251.974 486.187,217.339 506.933,210.796 527.68,226.853 548.426,228.72 569.173,226.126 589.92,236.753 610.666,232.586 631.413,203.849 \n",
       "  652.159,197.602 672.906,215.117 693.652,218.345 714.399,215.837 735.146,202.93 755.892,181.772 776.639,161.569 797.385,139.609 818.132,109.959 838.879,86.2547 \n",
       "  859.625,105.566 880.372,153.078 901.118,178.997 921.865,198.172 942.611,183.964 963.358,180.506 984.105,212.293 1004.85,185.137 1025.6,195.271 1046.34,193.603 \n",
       "  1067.09,176.413 1087.84,267.063 1108.58,263.187 1129.33,294.589 1150.08,307.354 1170.82,275.471 1191.57,320.787 1212.32,285.121 1233.06,234.26 1253.81,168.155 \n",
       "  1274.56,165.806 1295.3,244.191 1316.05,258.184 1336.8,253.438 1357.54,227.485 1378.29,148.896 1399.04,150.199 1419.78,150.592 1440.53,154.443 1461.28,141.343 \n",
       "  1482.02,132.125 1502.77,190.757 1523.52,186.804 1544.26,225.23 1565.01,255.47 1585.76,273.145 1606.5,337.922 1627.25,372.137 1647.99,417.061 1668.74,465.531 \n",
       "  1689.49,480.072 1710.23,522.581 1730.98,515.71 1751.73,547.132 1772.47,758.858 1793.22,763.302 1813.97,784.585 1834.71,794.947 1855.46,881.002 1876.21,933.266 \n",
       "  1896.95,966.67 1917.7,999.12 1938.45,1029.87 1959.19,1068.46 1979.94,1092.39 2000.69,1118.88 2021.43,1159.01 2042.18,1176.57 2062.93,1200.38 2083.67,1220.05 \n",
       "  2104.42,1244.19 2125.17,1271.02 2145.91,1290.64 2166.66,1312.15 2187.41,1331.25 2208.15,1346.33 2228.9,1360.62 2249.65,1372.01 2270.39,1379.76 2291.14,1386.61 \n",
       "  \n",
       "  \"/>\n",
       "</svg>\n"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Import 3 100-vectors: population, birth_rate, death_rate\n",
    "D = population_data();\n",
    "b = D[\"birth_rate\"];\n",
    "d = D[\"death_rate\"];\n",
    "A = [b'; diagonal(1 .- d[1:end-1]) zeros(length(d)-1)];\n",
    "x = D[\"population\"];\n",
    "for k = 1:10\n",
    "global x\n",
    "x = A*x;\n",
    "end;\n",
    "using Plots\n",
    "plot(x, legend=false, xlabel = \"Age\", ylabel = \"Population (millions)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Figure 9.2** Predicted age distribution in the US in 2020.\n",
    "\n",
    "\n",
    "Note the keyword `global` in the for-loop. Without this statement, the scope of the\n",
    "variable x created by the assignment `x = A*x` would be local to the for-loop, i.e.,\n",
    "this variable does not exist outside the loop and is different from the `x` outside the\n",
    "loop.\n",
    "\n",
    "**9.3 Epidemic dynamics**\n",
    "Let’s implement the simulation of the epidemic dynamics from VMLS §[9.3](https://web.stanford.edu/\\\\%7Eboyd/vmls/vmls.pdf#section.9.3). The plot is in figure 9.3."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/svg+xml": [
       "<?xml version=\"1.0\" encoding=\"utf-8\"?>\n",
       "<svg xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" width=\"600\" height=\"400\" viewBox=\"0 0 2400 1600\">\n",
       "<defs>\n",
       "  <clipPath id=\"clip5800\">\n",
       "    <rect x=\"0\" y=\"0\" width=\"2400\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip5800)\" d=\"\n",
       "M0 1600 L2400 1600 L2400 0 L0 0  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip5801\">\n",
       "    <rect x=\"480\" y=\"0\" width=\"1681\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip5800)\" d=\"\n",
       "M180.66 1425.62 L2352.76 1425.62 L2352.76 47.2441 L180.66 47.2441  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip5802\">\n",
       "    <rect x=\"180\" y=\"47\" width=\"2173\" height=\"1379\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  232.33,1425.62 232.33,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  722.557,1425.62 722.557,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1212.78,1425.62 1212.78,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1703.01,1425.62 1703.01,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  2193.24,1425.62 2193.24,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  180.66,1386.61 2352.76,1386.61 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  180.66,1061.52 2352.76,1061.52 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  180.66,736.431 2352.76,736.431 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  180.66,411.343 2352.76,411.343 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  180.66,86.2547 2352.76,86.2547 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  180.66,1425.62 2352.76,1425.62 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  180.66,1425.62 180.66,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  232.33,1425.62 232.33,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  722.557,1425.62 722.557,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1212.78,1425.62 1212.78,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1703.01,1425.62 1703.01,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  2193.24,1425.62 2193.24,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  180.66,1386.61 213.242,1386.61 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  180.66,1061.52 213.242,1061.52 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  180.66,736.431 213.242,736.431 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  180.66,411.343 213.242,411.343 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  180.66,86.2547 213.242,86.2547 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 232.33, 1479.62)\" x=\"232.33\" y=\"1479.62\">0</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 722.557, 1479.62)\" x=\"722.557\" y=\"1479.62\">50</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1212.78, 1479.62)\" x=\"1212.78\" y=\"1479.62\">100</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1703.01, 1479.62)\" x=\"1703.01\" y=\"1479.62\">150</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 2193.24, 1479.62)\" x=\"2193.24\" y=\"1479.62\">200</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 156.66, 1404.11)\" x=\"156.66\" y=\"1404.11\">0.00</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 156.66, 1079.02)\" x=\"156.66\" y=\"1079.02\">0.25</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 156.66, 753.931)\" x=\"156.66\" y=\"753.931\">0.50</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 156.66, 428.843)\" x=\"156.66\" y=\"428.843\">0.75</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 156.66, 103.755)\" x=\"156.66\" y=\"103.755\">1.00</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:66px; text-anchor:middle;\" transform=\"rotate(0, 1266.71, 1559.48)\" x=\"1266.71\" y=\"1559.48\">Time t</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#009af9; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.135,86.2547 251.939,151.272 261.744,210.438 271.548,264.566 281.353,314.336 291.157,360.324 300.962,403.012 310.766,442.807 320.571,480.052 330.375,515.04 \n",
       "  340.18,548.015 349.985,579.19 359.789,608.742 369.594,636.824 379.398,663.568 389.203,689.087 399.007,713.479 408.812,736.829 418.616,759.211 428.421,780.69 \n",
       "  438.225,801.325 448.03,821.165 457.834,840.256 467.639,858.639 477.443,876.35 487.248,893.422 497.053,909.887 506.857,925.771 516.662,941.101 526.466,955.899 \n",
       "  536.271,970.187 546.075,983.987 555.88,997.318 565.684,1010.2 575.489,1022.64 585.293,1034.67 595.098,1046.29 604.902,1057.52 614.707,1068.38 624.511,1078.88 \n",
       "  634.316,1089.02 644.12,1098.83 653.925,1108.32 663.73,1117.49 673.534,1126.36 683.339,1134.93 693.143,1143.22 702.948,1151.24 712.752,1158.99 722.557,1166.49 \n",
       "  732.361,1173.74 742.166,1180.75 751.97,1187.52 761.775,1194.08 771.579,1200.42 781.384,1206.55 791.188,1212.48 800.993,1218.21 810.797,1223.75 820.602,1229.12 \n",
       "  830.407,1234.3 840.211,1239.31 850.016,1244.16 859.82,1248.85 869.625,1253.39 879.429,1257.77 889.234,1262.01 899.038,1266.12 908.843,1270.08 918.647,1273.92 \n",
       "  928.452,1277.63 938.256,1281.21 948.061,1284.68 957.865,1288.04 967.67,1291.28 977.474,1294.42 987.279,1297.46 997.084,1300.39 1006.89,1303.23 1016.69,1305.97 \n",
       "  1026.5,1308.63 1036.3,1311.2 1046.11,1313.68 1055.91,1316.08 1065.72,1318.4 1075.52,1320.65 1085.32,1322.82 1095.13,1324.92 1104.93,1326.95 1114.74,1328.91 \n",
       "  1124.54,1330.81 1134.35,1332.65 1144.15,1334.42 1153.96,1336.14 1163.76,1337.8 1173.57,1339.41 1183.37,1340.96 1193.17,1342.47 1202.98,1343.92 1212.78,1345.32 \n",
       "  1222.59,1346.68 1232.39,1348 1242.2,1349.27 1252,1350.5 1261.81,1351.69 1271.61,1352.84 1281.41,1353.95 1291.22,1355.02 1301.02,1356.06 1310.83,1357.07 \n",
       "  1320.63,1358.04 1330.44,1358.98 1340.24,1359.89 1350.05,1360.77 1359.85,1361.62 1369.66,1362.44 1379.46,1363.24 1389.26,1364.01 1399.07,1364.75 1408.87,1365.47 \n",
       "  1418.68,1366.17 1428.48,1366.84 1438.29,1367.49 1448.09,1368.12 1457.9,1368.73 1467.7,1369.32 1477.51,1369.89 1487.31,1370.44 1497.11,1370.97 1506.92,1371.48 \n",
       "  1516.72,1371.98 1526.53,1372.46 1536.33,1372.93 1546.14,1373.38 1555.94,1373.81 1565.75,1374.24 1575.55,1374.64 1585.36,1375.04 1595.16,1375.42 1604.96,1375.79 \n",
       "  1614.77,1376.14 1624.57,1376.49 1634.38,1376.82 1644.18,1377.14 1653.99,1377.45 1663.79,1377.75 1673.6,1378.05 1683.4,1378.33 1693.21,1378.6 1703.01,1378.86 \n",
       "  1712.81,1379.12 1722.62,1379.37 1732.42,1379.6 1742.23,1379.83 1752.03,1380.06 1761.84,1380.27 1771.64,1380.48 1781.45,1380.68 1791.25,1380.88 1801.05,1381.07 \n",
       "  1810.86,1381.25 1820.66,1381.43 1830.47,1381.6 1840.27,1381.76 1850.08,1381.92 1859.88,1382.08 1869.69,1382.22 1879.49,1382.37 1889.3,1382.51 1899.1,1382.64 \n",
       "  1908.9,1382.77 1918.71,1382.9 1928.51,1383.02 1938.32,1383.14 1948.12,1383.25 1957.93,1383.36 1967.73,1383.47 1977.54,1383.57 1987.34,1383.67 1997.15,1383.77 \n",
       "  2006.95,1383.86 2016.75,1383.95 2026.56,1384.04 2036.36,1384.13 2046.17,1384.21 2055.97,1384.29 2065.78,1384.36 2075.58,1384.44 2085.39,1384.51 2095.19,1384.58 \n",
       "  2105,1384.64 2114.8,1384.71 2124.6,1384.77 2134.41,1384.83 2144.21,1384.89 2154.02,1384.95 2163.82,1385 2173.63,1385.05 2183.43,1385.11 2193.24,1385.16 \n",
       "  2203.04,1385.2 2212.85,1385.25 2222.65,1385.29 2232.45,1385.34 2242.26,1385.38 2252.06,1385.42 2261.87,1385.46 2271.67,1385.5 2281.48,1385.53 2291.28,1385.57 \n",
       "  \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#e26f46; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.135,1386.61 251.939,1321.59 261.744,1269.58 271.548,1228.32 281.353,1195.96 291.157,1170.95 300.962,1151.98 310.766,1138 320.571,1128.1 330.375,1121.55 \n",
       "  340.18,1117.73 349.985,1116.13 359.789,1116.33 369.594,1117.98 379.398,1120.78 389.203,1124.51 399.007,1128.94 408.812,1133.94 418.616,1139.35 428.421,1145.07 \n",
       "  438.225,1151 448.03,1157.08 457.834,1163.24 467.639,1169.43 477.443,1175.6 487.248,1181.74 497.053,1187.81 506.857,1193.8 516.662,1199.68 526.466,1205.44 \n",
       "  536.271,1211.08 546.075,1216.59 555.88,1221.96 565.684,1227.19 575.489,1232.28 585.293,1237.23 595.098,1242.04 604.902,1246.71 614.707,1251.24 624.511,1255.64 \n",
       "  634.316,1259.89 644.12,1264.02 653.925,1268.02 663.73,1271.9 673.534,1275.65 683.339,1279.28 693.143,1282.79 702.948,1286.2 712.752,1289.49 722.557,1292.68 \n",
       "  732.361,1295.76 742.166,1298.74 751.97,1301.63 761.775,1304.42 771.579,1307.12 781.384,1309.74 791.188,1312.27 800.993,1314.71 810.797,1317.07 820.602,1319.36 \n",
       "  830.407,1321.57 840.211,1323.71 850.016,1325.78 859.82,1327.78 869.625,1329.72 879.429,1331.59 889.234,1333.4 899.038,1335.15 908.843,1336.85 918.647,1338.49 \n",
       "  928.452,1340.07 938.256,1341.6 948.061,1343.08 957.865,1344.52 967.67,1345.9 977.474,1347.24 987.279,1348.54 997.084,1349.79 1006.89,1351 1016.69,1352.17 \n",
       "  1026.5,1353.31 1036.3,1354.4 1046.11,1355.46 1055.91,1356.49 1065.72,1357.48 1075.52,1358.44 1085.32,1359.37 1095.13,1360.26 1104.93,1361.13 1114.74,1361.97 \n",
       "  1124.54,1362.78 1134.35,1363.56 1144.15,1364.32 1153.96,1365.06 1163.76,1365.77 1173.57,1366.45 1183.37,1367.12 1193.17,1367.76 1202.98,1368.38 1212.78,1368.98 \n",
       "  1222.59,1369.56 1232.39,1370.12 1242.2,1370.66 1252,1371.19 1261.81,1371.69 1271.61,1372.19 1281.41,1372.66 1291.22,1373.12 1301.02,1373.56 1310.83,1373.99 \n",
       "  1320.63,1374.41 1330.44,1374.81 1340.24,1375.2 1350.05,1375.57 1359.85,1375.94 1369.66,1376.29 1379.46,1376.63 1389.26,1376.96 1399.07,1377.27 1408.87,1377.58 \n",
       "  1418.68,1377.88 1428.48,1378.17 1438.29,1378.44 1448.09,1378.71 1457.9,1378.97 1467.7,1379.22 1477.51,1379.47 1487.31,1379.7 1497.11,1379.93 1506.92,1380.15 \n",
       "  1516.72,1380.36 1526.53,1380.57 1536.33,1380.77 1546.14,1380.96 1555.94,1381.14 1565.75,1381.32 1575.55,1381.5 1585.36,1381.67 1595.16,1381.83 1604.96,1381.99 \n",
       "  1614.77,1382.14 1624.57,1382.29 1634.38,1382.43 1644.18,1382.57 1653.99,1382.7 1663.79,1382.83 1673.6,1382.95 1683.4,1383.07 1693.21,1383.19 1703.01,1383.3 \n",
       "  1712.81,1383.41 1722.62,1383.52 1732.42,1383.62 1742.23,1383.72 1752.03,1383.81 1761.84,1383.9 1771.64,1383.99 1781.45,1384.08 1791.25,1384.16 1801.05,1384.24 \n",
       "  1810.86,1384.32 1820.66,1384.39 1830.47,1384.47 1840.27,1384.54 1850.08,1384.61 1859.88,1384.67 1869.69,1384.74 1879.49,1384.8 1889.3,1384.86 1899.1,1384.91 \n",
       "  1908.9,1384.97 1918.71,1385.02 1928.51,1385.08 1938.32,1385.13 1948.12,1385.18 1957.93,1385.22 1967.73,1385.27 1977.54,1385.31 1987.34,1385.36 1997.15,1385.4 \n",
       "  2006.95,1385.44 2016.75,1385.47 2026.56,1385.51 2036.36,1385.55 2046.17,1385.58 2055.97,1385.62 2065.78,1385.65 2075.58,1385.68 2085.39,1385.71 2095.19,1385.74 \n",
       "  2105,1385.77 2114.8,1385.8 2124.6,1385.82 2134.41,1385.85 2144.21,1385.87 2154.02,1385.9 2163.82,1385.92 2173.63,1385.94 2183.43,1385.97 2193.24,1385.99 \n",
       "  2203.04,1386.01 2212.85,1386.03 2222.65,1386.05 2232.45,1386.07 2242.26,1386.08 2252.06,1386.1 2261.87,1386.12 2271.67,1386.13 2281.48,1386.15 2291.28,1386.16 \n",
       "  \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#3da44d; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.135,1386.61 251.939,1386.61 261.744,1380.11 271.548,1368.4 281.353,1352.57 291.157,1333.51 300.962,1311.94 310.766,1288.48 320.571,1263.62 330.375,1237.77 \n",
       "  340.18,1211.26 349.985,1184.37 359.789,1157.33 369.594,1130.3 379.398,1103.44 389.203,1076.85 399.007,1050.64 408.812,1024.88 418.616,999.61 428.421,974.884 \n",
       "  438.225,950.73 448.03,927.17 457.834,904.217 467.639,881.88 477.443,860.162 487.248,839.061 497.053,818.575 506.857,798.695 516.662,779.414 526.466,760.721 \n",
       "  536.271,742.604 546.075,725.051 555.88,708.049 565.684,691.584 575.489,675.643 585.293,660.21 595.098,645.273 604.902,630.816 614.707,616.827 624.511,603.29 \n",
       "  634.316,590.193 644.12,577.522 653.925,565.263 663.73,553.405 673.534,541.933 683.339,530.837 693.143,520.104 702.948,509.723 712.752,499.682 722.557,489.97 \n",
       "  732.361,480.577 742.166,471.492 751.97,462.706 761.775,454.208 771.579,445.989 781.384,438.041 791.188,430.354 800.993,422.92 810.797,415.73 820.602,408.777 \n",
       "  830.407,402.052 840.211,395.549 850.016,389.259 859.82,383.177 869.625,377.294 879.429,371.606 889.234,366.104 899.038,360.784 908.843,355.638 918.647,350.662 \n",
       "  928.452,345.85 938.256,341.196 948.061,336.695 957.865,332.343 967.67,328.134 977.474,324.063 987.279,320.126 997.084,316.319 1006.89,312.637 1016.69,309.076 \n",
       "  1026.5,305.633 1036.3,302.303 1046.11,299.082 1055.91,295.968 1065.72,292.956 1075.52,290.043 1085.32,287.226 1095.13,284.502 1104.93,281.867 1114.74,279.32 \n",
       "  1124.54,276.856 1134.35,274.473 1144.15,272.168 1153.96,269.94 1163.76,267.785 1173.57,265.7 1183.37,263.685 1193.17,261.736 1202.98,259.85 1212.78,258.027 \n",
       "  1222.59,256.264 1232.39,254.559 1242.2,252.91 1252,251.316 1261.81,249.774 1271.61,248.282 1281.41,246.84 1291.22,245.445 1301.02,244.096 1310.83,242.792 \n",
       "  1320.63,241.53 1330.44,240.31 1340.24,239.131 1350.05,237.99 1359.85,236.886 1369.66,235.819 1379.46,234.787 1389.26,233.789 1399.07,232.824 1408.87,231.89 \n",
       "  1418.68,230.988 1428.48,230.115 1438.29,229.271 1448.09,228.454 1457.9,227.665 1467.7,226.901 1477.51,226.163 1487.31,225.449 1497.11,224.758 1506.92,224.09 \n",
       "  1516.72,223.444 1526.53,222.82 1536.33,222.215 1546.14,221.631 1555.94,221.066 1565.75,220.52 1575.55,219.992 1585.36,219.481 1595.16,218.986 1604.96,218.509 \n",
       "  1614.77,218.046 1624.57,217.599 1634.38,217.167 1644.18,216.749 1653.99,216.345 1663.79,215.954 1673.6,215.576 1683.4,215.21 1693.21,214.857 1703.01,214.515 \n",
       "  1712.81,214.184 1722.62,213.864 1732.42,213.555 1742.23,213.256 1752.03,212.967 1761.84,212.687 1771.64,212.416 1781.45,212.155 1791.25,211.902 1801.05,211.657 \n",
       "  1810.86,211.42 1820.66,211.192 1830.47,210.97 1840.27,210.756 1850.08,210.549 1859.88,210.349 1869.69,210.155 1879.49,209.968 1889.3,209.787 1899.1,209.612 \n",
       "  1908.9,209.443 1918.71,209.279 1928.51,209.121 1938.32,208.968 1948.12,208.82 1957.93,208.676 1967.73,208.538 1977.54,208.404 1987.34,208.274 1997.15,208.149 \n",
       "  2006.95,208.028 2016.75,207.911 2026.56,207.797 2036.36,207.688 2046.17,207.582 2055.97,207.479 2065.78,207.38 2075.58,207.284 2085.39,207.192 2095.19,207.102 \n",
       "  2105,207.015 2114.8,206.932 2124.6,206.851 2134.41,206.772 2144.21,206.696 2154.02,206.623 2163.82,206.552 2173.63,206.483 2183.43,206.417 2193.24,206.353 \n",
       "  2203.04,206.291 2212.85,206.231 2222.65,206.173 2232.45,206.117 2242.26,206.063 2252.06,206.01 2261.87,205.959 2271.67,205.91 2281.48,205.863 2291.28,205.817 \n",
       "  \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5802)\" style=\"stroke:#c271d2; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.135,1386.61 251.939,1386.61 261.744,1385.96 271.548,1384.79 281.353,1383.2 291.157,1381.3 300.962,1379.14 310.766,1376.8 320.571,1374.31 330.375,1371.72 \n",
       "  340.18,1369.07 349.985,1366.38 359.789,1363.68 369.594,1360.98 379.398,1358.29 389.203,1355.63 399.007,1353.01 408.812,1350.43 418.616,1347.91 428.421,1345.44 \n",
       "  438.225,1343.02 448.03,1340.66 457.834,1338.37 467.639,1336.14 477.443,1333.96 487.248,1331.85 497.053,1329.8 506.857,1327.82 516.662,1325.89 526.466,1324.02 \n",
       "  536.271,1322.21 546.075,1320.45 555.88,1318.75 565.684,1317.11 575.489,1315.51 585.293,1313.97 595.098,1312.47 604.902,1311.03 614.707,1309.63 624.511,1308.28 \n",
       "  634.316,1306.97 644.12,1305.7 653.925,1304.47 663.73,1303.29 673.534,1302.14 683.339,1301.03 693.143,1299.96 702.948,1298.92 712.752,1297.92 722.557,1296.94 \n",
       "  732.361,1296 742.166,1295.1 751.97,1294.22 761.775,1293.37 771.579,1292.55 781.384,1291.75 791.188,1290.98 800.993,1290.24 810.797,1289.52 820.602,1288.82 \n",
       "  830.407,1288.15 840.211,1287.5 850.016,1286.87 859.82,1286.26 869.625,1285.68 879.429,1285.11 889.234,1284.56 899.038,1284.03 908.843,1283.51 918.647,1283.01 \n",
       "  928.452,1282.53 938.256,1282.07 948.061,1281.62 957.865,1281.18 967.67,1280.76 977.474,1280.35 987.279,1279.96 997.084,1279.58 1006.89,1279.21 1016.69,1278.85 \n",
       "  1026.5,1278.51 1036.3,1278.18 1046.11,1277.86 1055.91,1277.54 1065.72,1277.24 1075.52,1276.95 1085.32,1276.67 1095.13,1276.4 1104.93,1276.13 1114.74,1275.88 \n",
       "  1124.54,1275.63 1134.35,1275.39 1144.15,1275.16 1153.96,1274.94 1163.76,1274.73 1173.57,1274.52 1183.37,1274.32 1193.17,1274.12 1202.98,1273.93 1212.78,1273.75 \n",
       "  1222.59,1273.57 1232.39,1273.4 1242.2,1273.24 1252,1273.08 1261.81,1272.92 1271.61,1272.78 1281.41,1272.63 1291.22,1272.49 1301.02,1272.36 1310.83,1272.23 \n",
       "  1320.63,1272.1 1330.44,1271.98 1340.24,1271.86 1350.05,1271.75 1359.85,1271.64 1369.66,1271.53 1379.46,1271.43 1389.26,1271.33 1399.07,1271.23 1408.87,1271.14 \n",
       "  1418.68,1271.05 1428.48,1270.96 1438.29,1270.87 1448.09,1270.79 1457.9,1270.71 1467.7,1270.64 1477.51,1270.56 1487.31,1270.49 1497.11,1270.42 1506.92,1270.36 \n",
       "  1516.72,1270.29 1526.53,1270.23 1536.33,1270.17 1546.14,1270.11 1555.94,1270.05 1565.75,1270 1575.55,1269.95 1585.36,1269.9 1595.16,1269.85 1604.96,1269.8 \n",
       "  1614.77,1269.75 1624.57,1269.71 1634.38,1269.66 1644.18,1269.62 1653.99,1269.58 1663.79,1269.54 1673.6,1269.5 1683.4,1269.47 1693.21,1269.43 1703.01,1269.4 \n",
       "  1712.81,1269.37 1722.62,1269.33 1732.42,1269.3 1742.23,1269.27 1752.03,1269.24 1761.84,1269.22 1771.64,1269.19 1781.45,1269.16 1791.25,1269.14 1801.05,1269.11 \n",
       "  1810.86,1269.09 1820.66,1269.07 1830.47,1269.04 1840.27,1269.02 1850.08,1269 1859.88,1268.98 1869.69,1268.96 1879.49,1268.94 1889.3,1268.93 1899.1,1268.91 \n",
       "  1908.9,1268.89 1918.71,1268.88 1928.51,1268.86 1938.32,1268.84 1948.12,1268.83 1957.93,1268.81 1967.73,1268.8 1977.54,1268.79 1987.34,1268.77 1997.15,1268.76 \n",
       "  2006.95,1268.75 2016.75,1268.74 2026.56,1268.73 2036.36,1268.72 2046.17,1268.71 2055.97,1268.7 2065.78,1268.69 2075.58,1268.68 2085.39,1268.67 2095.19,1268.66 \n",
       "  2105,1268.65 2114.8,1268.64 2124.6,1268.63 2134.41,1268.62 2144.21,1268.62 2154.02,1268.61 2163.82,1268.6 2173.63,1268.6 2183.43,1268.59 2193.24,1268.58 \n",
       "  2203.04,1268.58 2212.85,1268.57 2222.65,1268.56 2232.45,1268.56 2242.26,1268.55 2252.06,1268.55 2261.87,1268.54 2271.67,1268.54 2281.48,1268.53 2291.28,1268.53 \n",
       "  \n",
       "  \"/>\n",
       "<path clip-path=\"url(#clip5800)\" d=\"\n",
       "M1791.95 433.164 L2280.76 433.164 L2280.76 130.764 L1791.95 130.764  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1791.95,433.164 2280.76,433.164 2280.76,130.764 1791.95,130.764 1791.95,433.164 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#009af9; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1815.95,191.244 1959.95,191.244 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:start;\" transform=\"rotate(0, 1983.95, 208.744)\" x=\"1983.95\" y=\"208.744\">Susceptible</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#e26f46; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1815.95,251.724 1959.95,251.724 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:start;\" transform=\"rotate(0, 1983.95, 269.224)\" x=\"1983.95\" y=\"269.224\">Infected</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#3da44d; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1815.95,312.204 1959.95,312.204 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:start;\" transform=\"rotate(0, 1983.95, 329.704)\" x=\"1983.95\" y=\"329.704\">Recovered</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip5800)\" style=\"stroke:#c271d2; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1815.95,372.684 1959.95,372.684 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip5800)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:start;\" transform=\"rotate(0, 1983.95, 390.184)\" x=\"1983.95\" y=\"390.184\">Deceased</text>\n",
       "</g>\n",
       "</svg>\n"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "T = 210;\n",
    "A = [ 0.95 0.04 0 0 ; 0.05 0.85 0 0 ; 0 0.10 1 0 ; 0 0.01 0 1 ];\n",
    "x_1 = [1,0,0,0];\n",
    "state_traj = [x_1 zeros(4,T-1) ]; # State trajectory\n",
    "for t=1:T-1 # Dynamics recursion\n",
    "state_traj[:,t+1] = A*state_traj[:,t];\n",
    "end\n",
    "using Plots\n",
    "plot(1:T, state_traj', xlabel = \"Time t\", label = [\"Susceptible\", \"Infected\", \"Recovered\", \"Deceased\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Figure 9.3** Simulation of epidemic dynamics."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**9.4 Motion of a mass**\n",
    "Let’s simulate the discretized model of the motion of a mass in §[9.4](https://web.stanford.edu/~boyd/vmls/vmls.pdf#section.9.4) of VMLS. See figure 9.4."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/svg+xml": [
       "<?xml version=\"1.0\" encoding=\"utf-8\"?>\n",
       "<svg xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" width=\"600\" height=\"400\" viewBox=\"0 0 2400 1600\">\n",
       "<defs>\n",
       "  <clipPath id=\"clip6600\">\n",
       "    <rect x=\"0\" y=\"0\" width=\"2400\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip6600)\" d=\"\n",
       "M0 1600 L2400 1600 L2400 0 L0 0  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip6601\">\n",
       "    <rect x=\"480\" y=\"0\" width=\"1681\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip6600)\" d=\"\n",
       "M242.516 1425.62 L2352.76 1425.62 L2352.76 47.2441 L242.516 47.2441  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip6602\">\n",
       "    <rect x=\"242\" y=\"47\" width=\"2111\" height=\"1379\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  298.917,1425.62 298.917,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  631.269,1425.62 631.269,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  963.622,1425.62 963.622,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1295.97,1425.62 1295.97,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1628.33,1425.62 1628.33,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1960.68,1425.62 1960.68,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  2293.03,1425.62 2293.03,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  242.516,1250.28 2352.76,1250.28 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  242.516,879.175 2352.76,879.175 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  242.516,508.075 2352.76,508.075 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  242.516,136.975 2352.76,136.975 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.516,1425.62 2352.76,1425.62 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.516,1425.62 242.516,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  298.917,1425.62 298.917,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  631.269,1425.62 631.269,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  963.622,1425.62 963.622,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1295.97,1425.62 1295.97,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1628.33,1425.62 1628.33,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1960.68,1425.62 1960.68,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  2293.03,1425.62 2293.03,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.516,1250.28 274.17,1250.28 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.516,879.175 274.17,879.175 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.516,508.075 274.17,508.075 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip6600)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  242.516,136.975 274.17,136.975 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 298.917, 1479.62)\" x=\"298.917\" y=\"1479.62\">0</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 631.269, 1479.62)\" x=\"631.269\" y=\"1479.62\">100</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 963.622, 1479.62)\" x=\"963.622\" y=\"1479.62\">200</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1295.97, 1479.62)\" x=\"1295.97\" y=\"1479.62\">300</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1628.33, 1479.62)\" x=\"1628.33\" y=\"1479.62\">400</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1960.68, 1479.62)\" x=\"1960.68\" y=\"1479.62\">500</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 2293.03, 1479.62)\" x=\"2293.03\" y=\"1479.62\">600</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 218.516, 1267.78)\" x=\"218.516\" y=\"1267.78\">0.00</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 218.516, 896.675)\" x=\"218.516\" y=\"896.675\">0.05</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 218.516, 525.575)\" x=\"218.516\" y=\"525.575\">0.10</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 218.516, 154.475)\" x=\"218.516\" y=\"154.475\">0.15</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:66px; text-anchor:middle;\" transform=\"rotate(0, 1297.64, 1559.48)\" x=\"1297.64\" y=\"1559.48\">k</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip6600)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:66px; text-anchor:middle;\" transform=\"rotate(-90, 89.2861, 736.431)\" x=\"89.2861\" y=\"736.431\">Position</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip6602)\" style=\"stroke:#009af9; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  302.24,1250.28 305.564,1250.28 308.887,1250.28 312.211,1250.28 315.534,1250.28 318.858,1250.28 322.181,1250.28 325.505,1250.28 328.828,1250.28 332.152,1250.28 \n",
       "  335.475,1250.28 338.799,1250.28 342.122,1250.28 345.446,1250.28 348.77,1250.28 352.093,1250.28 355.417,1250.28 358.74,1250.28 362.064,1250.28 365.387,1250.28 \n",
       "  368.711,1250.28 372.034,1250.28 375.358,1250.28 378.681,1250.28 382.005,1250.28 385.328,1250.28 388.652,1250.28 391.975,1250.28 395.299,1250.28 398.622,1250.28 \n",
       "  401.946,1250.28 405.269,1250.28 408.593,1250.28 411.917,1250.28 415.24,1250.28 418.564,1250.28 421.887,1250.28 425.211,1250.28 428.534,1250.28 431.858,1250.28 \n",
       "  435.181,1250.28 438.505,1250.28 441.828,1250.28 445.152,1250.28 448.475,1250.28 451.799,1250.28 455.122,1250.28 458.446,1250.28 461.769,1250.28 465.093,1250.28 \n",
       "  468.416,1250.28 471.74,1249.53 475.064,1248.06 478.387,1245.85 481.711,1242.93 485.034,1239.29 488.358,1234.95 491.681,1229.9 495.005,1224.17 498.328,1217.75 \n",
       "  501.652,1210.66 504.975,1202.89 508.299,1194.45 511.622,1185.36 514.946,1175.62 518.269,1165.24 521.593,1154.21 524.916,1142.55 528.24,1130.27 531.563,1117.37 \n",
       "  534.887,1103.86 538.211,1089.73 541.534,1075.01 544.858,1059.69 548.181,1043.79 551.505,1027.3 554.828,1010.23 558.152,992.589 561.475,974.384 564.799,955.619 \n",
       "  568.122,936.3 571.446,916.431 574.769,896.019 578.093,875.069 581.416,853.587 584.74,831.577 588.063,809.044 591.387,785.995 594.71,762.434 598.034,738.367 \n",
       "  601.358,713.798 604.681,688.733 608.005,663.176 611.328,637.132 614.652,610.606 617.975,583.604 621.299,556.13 624.622,528.188 627.946,499.783 631.269,470.92 \n",
       "  634.593,441.604 637.916,413.545 641.24,386.732 644.563,361.152 647.887,336.793 651.21,313.642 654.534,291.687 657.857,270.917 661.181,251.32 664.504,232.883 \n",
       "  667.828,215.595 671.152,199.446 674.475,184.422 677.799,170.514 681.122,157.709 684.446,145.998 687.769,135.369 691.093,125.81 694.416,117.312 697.74,109.864 \n",
       "  701.063,103.456 704.387,98.0758 707.71,93.7147 711.034,90.3621 714.357,88.0078 717.681,86.642 721.004,86.2547 724.328,86.8361 727.651,88.3766 730.975,90.8665 \n",
       "  734.299,94.2964 737.622,98.6568 740.946,103.939 744.269,110.132 747.593,117.229 750.916,125.22 754.24,134.095 757.563,143.847 760.887,154.466 764.21,165.943 \n",
       "  767.534,178.271 770.857,190.475 774.181,202.558 777.504,214.519 780.828,226.361 784.151,238.085 787.475,249.691 790.798,261.182 794.122,272.557 797.446,283.818 \n",
       "  800.769,294.967 804.093,306.005 807.416,316.932 810.74,327.75 814.063,338.46 817.387,349.062 820.71,359.559 824.034,369.95 827.357,380.238 830.681,390.423 \n",
       "  834.004,400.506 837.328,410.488 840.651,420.37 843.975,430.153 847.298,439.839 850.622,449.428 853.945,458.921 857.269,468.319 860.593,477.623 863.916,486.834 \n",
       "  867.24,495.952 870.563,504.98 873.887,513.917 877.21,522.765 880.534,531.525 883.857,540.197 887.181,548.782 890.504,557.281 893.828,565.696 897.151,574.026 \n",
       "  900.475,582.273 903.798,590.437 907.122,598.52 910.445,606.522 913.769,614.444 917.092,622.287 920.416,630.051 923.74,637.738 927.063,645.347 930.387,652.881 \n",
       "  933.71,660.339 937.034,667.723 940.357,675.033 943.681,682.27 947.004,689.434 950.328,696.527 953.651,703.549 956.975,710.501 960.298,717.383 963.622,724.196 \n",
       "  966.945,730.941 970.269,737.619 973.592,744.23 976.916,750.775 980.239,757.254 983.563,763.669 986.887,770.019 990.21,776.306 993.534,782.531 996.857,788.692 \n",
       "  1000.18,794.793 1003.5,800.832 1006.83,806.811 1010.15,812.73 1013.47,818.59 1016.8,824.391 1020.12,830.134 1023.45,835.82 1026.77,841.449 1030.09,847.022 \n",
       "  1033.42,852.539 1036.74,858 1040.06,863.407 1043.39,868.761 1046.71,874.06 1050.03,879.307 1053.36,884.501 1056.68,889.643 1060,894.734 1063.33,899.773 \n",
       "  1066.65,904.763 1069.97,909.702 1073.3,914.593 1076.62,919.434 1079.95,924.227 1083.27,928.972 1086.59,933.669 1089.92,938.319 1093.24,942.923 1096.56,947.481 \n",
       "  1099.89,951.994 1103.21,956.461 1106.53,960.883 1109.86,965.262 1113.18,969.596 1116.5,973.888 1119.83,978.136 1123.15,982.342 1126.47,986.505 1129.8,990.627 \n",
       "  1133.12,994.708 1136.45,998.748 1139.77,1002.75 1143.09,1006.71 1146.42,1010.63 1149.74,1014.51 1153.06,1018.35 1156.39,1022.15 1159.71,1025.92 1163.03,1029.65 \n",
       "  1166.36,1033.34 1169.68,1036.99 1173,1040.61 1176.33,1044.19 1179.65,1047.74 1182.97,1051.25 1186.3,1054.72 1189.62,1058.16 1192.95,1061.57 1196.27,1064.94 \n",
       "  1199.59,1068.28 1202.92,1071.58 1206.24,1074.85 1209.56,1078.09 1212.89,1081.3 1216.21,1084.47 1219.53,1087.61 1222.86,1090.72 1226.18,1093.8 1229.5,1096.85 \n",
       "  1232.83,1099.87 1236.15,1102.86 1239.47,1105.82 1242.8,1108.75 1246.12,1111.65 1249.45,1114.52 1252.77,1117.36 1256.09,1120.17 1259.42,1122.96 1262.74,1125.72 \n",
       "  1266.06,1128.45 1269.39,1131.15 1272.71,1133.82 1276.03,1136.47 1279.36,1139.1 1282.68,1141.69 1286,1144.26 1289.33,1146.81 1292.65,1149.33 1295.97,1151.82 \n",
       "  1299.3,1154.29 1302.62,1156.73 1305.94,1159.15 1309.27,1161.55 1312.59,1163.92 1315.92,1166.27 1319.24,1168.59 1322.56,1170.89 1325.89,1173.17 1329.21,1175.43 \n",
       "  1332.53,1177.66 1335.86,1179.87 1339.18,1182.06 1342.5,1184.23 1345.83,1186.37 1349.15,1188.49 1352.47,1190.6 1355.8,1192.68 1359.12,1194.74 1362.44,1196.78 \n",
       "  1365.77,1198.8 1369.09,1200.8 1372.42,1202.78 1375.74,1204.74 1379.06,1206.67 1382.39,1208.6 1385.71,1210.5 1389.03,1212.38 1392.36,1214.24 1395.68,1216.09 \n",
       "  1399,1217.91 1402.33,1219.72 1405.65,1221.51 1408.97,1223.28 1412.3,1225.04 1415.62,1226.77 1418.94,1228.49 1422.27,1230.2 1425.59,1231.88 1428.92,1233.55 \n",
       "  1432.24,1235.2 1435.56,1236.84 1438.89,1238.46 1442.21,1240.06 1445.53,1241.64 1448.86,1243.22 1452.18,1244.77 1455.5,1246.31 1458.83,1247.83 1462.15,1249.34 \n",
       "  1465.47,1250.84 1468.8,1252.31 1472.12,1253.78 1475.44,1255.23 1478.77,1256.66 1482.09,1258.08 1485.42,1259.49 1488.74,1260.88 1492.06,1262.26 1495.39,1263.63 \n",
       "  1498.71,1264.98 1502.03,1266.31 1505.36,1267.64 1508.68,1268.95 1512,1270.25 1515.33,1271.53 1518.65,1272.8 1521.97,1274.06 1525.3,1275.31 1528.62,1276.54 \n",
       "  1531.94,1277.76 1535.27,1278.97 1538.59,1280.17 1541.92,1281.36 1545.24,1282.53 1548.56,1283.69 1551.89,1284.84 1555.21,1285.98 1558.53,1287.11 1561.86,1288.22 \n",
       "  1565.18,1289.33 1568.5,1290.42 1571.83,1291.51 1575.15,1292.58 1578.47,1293.64 1581.8,1294.69 1585.12,1295.73 1588.44,1296.76 1591.77,1297.78 1595.09,1298.79 \n",
       "  1598.42,1299.79 1601.74,1300.78 1605.06,1301.76 1608.39,1302.73 1611.71,1303.69 1615.03,1304.64 1618.36,1305.58 1621.68,1306.51 1625,1307.43 1628.33,1308.34 \n",
       "  1631.65,1309.25 1634.97,1310.14 1638.3,1311.03 1641.62,1311.91 1644.94,1312.77 1648.27,1313.63 1651.59,1314.48 1654.92,1315.33 1658.24,1316.16 1661.56,1316.99 \n",
       "  1664.89,1317.8 1668.21,1318.61 1671.53,1319.41 1674.86,1320.21 1678.18,1320.99 1681.5,1321.77 1684.83,1322.54 1688.15,1323.3 1691.47,1324.05 1694.8,1324.8 \n",
       "  1698.12,1325.54 1701.44,1326.27 1704.77,1327 1708.09,1327.71 1711.42,1328.42 1714.74,1329.13 1718.06,1329.82 1721.39,1330.51 1724.71,1331.19 1728.03,1331.87 \n",
       "  1731.36,1332.54 1734.68,1333.2 1738,1333.85 1741.33,1334.5 1744.65,1335.14 1747.97,1335.78 1751.3,1336.41 1754.62,1337.03 1757.94,1337.65 1761.27,1338.26 \n",
       "  1764.59,1338.86 1767.92,1339.46 1771.24,1340.06 1774.56,1340.64 1777.89,1341.22 1781.21,1341.8 1784.53,1342.37 1787.86,1342.93 1791.18,1343.49 1794.5,1344.04 \n",
       "  1797.83,1344.59 1801.15,1345.13 1804.47,1345.66 1807.8,1346.19 1811.12,1346.72 1814.44,1347.24 1817.77,1347.75 1821.09,1348.26 1824.42,1348.77 1827.74,1349.27 \n",
       "  1831.06,1349.76 1834.39,1350.25 1837.71,1350.74 1841.03,1351.22 1844.36,1351.69 1847.68,1352.16 1851,1352.63 1854.33,1353.09 1857.65,1353.54 1860.97,1354 \n",
       "  1864.3,1354.44 1867.62,1354.89 1870.94,1355.32 1874.27,1355.76 1877.59,1356.19 1880.91,1356.61 1884.24,1357.03 1887.56,1357.45 1890.89,1357.86 1894.21,1358.27 \n",
       "  1897.53,1358.68 1900.86,1359.08 1904.18,1359.47 1907.5,1359.87 1910.83,1360.25 1914.15,1360.64 1917.47,1361.02 1920.8,1361.4 1924.12,1361.77 1927.44,1362.14 \n",
       "  1930.77,1362.51 1934.09,1362.87 1937.41,1363.23 1940.74,1363.58 1944.06,1363.93 1947.39,1364.28 1950.71,1364.62 1954.03,1364.97 1957.36,1365.3 1960.68,1365.64 \n",
       "  1964,1365.97 1967.33,1366.3 1970.65,1366.62 1973.97,1366.94 1977.3,1367.26 1980.62,1367.57 1983.94,1367.88 1987.27,1368.19 1990.59,1368.5 1993.91,1368.8 \n",
       "  1997.24,1369.1 2000.56,1369.4 2003.89,1369.69 2007.21,1369.98 2010.53,1370.27 2013.86,1370.55 2017.18,1370.83 2020.5,1371.11 2023.83,1371.39 2027.15,1371.66 \n",
       "  2030.47,1371.93 2033.8,1372.2 2037.12,1372.46 2040.44,1372.73 2043.77,1372.99 2047.09,1373.24 2050.41,1373.5 2053.74,1373.75 2057.06,1374 2060.39,1374.25 \n",
       "  2063.71,1374.49 2067.03,1374.73 2070.36,1374.97 2073.68,1375.21 2077,1375.45 2080.33,1375.68 2083.65,1375.91 2086.97,1376.14 2090.3,1376.36 2093.62,1376.59 \n",
       "  2096.94,1376.81 2100.27,1377.03 2103.59,1377.24 2106.91,1377.46 2110.24,1377.67 2113.56,1377.88 2116.89,1378.09 2120.21,1378.3 2123.53,1378.5 2126.86,1378.7 \n",
       "  2130.18,1378.9 2133.5,1379.1 2136.83,1379.3 2140.15,1379.49 2143.47,1379.68 2146.8,1379.87 2150.12,1380.06 2153.44,1380.25 2156.77,1380.43 2160.09,1380.62 \n",
       "  2163.41,1380.8 2166.74,1380.98 2170.06,1381.15 2173.39,1381.33 2176.71,1381.5 2180.03,1381.68 2183.36,1381.85 2186.68,1382.01 2190,1382.18 2193.33,1382.35 \n",
       "  2196.65,1382.51 2199.97,1382.67 2203.3,1382.83 2206.62,1382.99 2209.94,1383.15 2213.27,1383.31 2216.59,1383.46 2219.91,1383.61 2223.24,1383.76 2226.56,1383.91 \n",
       "  2229.89,1384.06 2233.21,1384.21 2236.53,1384.35 2239.86,1384.5 2243.18,1384.64 2246.5,1384.78 2249.83,1384.92 2253.15,1385.06 2256.47,1385.19 2259.8,1385.33 \n",
       "  2263.12,1385.46 2266.44,1385.59 2269.77,1385.73 2273.09,1385.86 2276.41,1385.98 2279.74,1386.11 2283.06,1386.24 2286.39,1386.36 2289.71,1386.49 2293.03,1386.61 \n",
       "  \n",
       "  \"/>\n",
       "</svg>\n"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "h = 0.01; m = 1; eta = 1;\n",
    "A = [ 1 h ; 0 1-h*eta/m ];\n",
    "B = [ 0 ; h/m ];\n",
    "x1 = [0,0];\n",
    "K = 600; # simulate for K*h = 6 seconds\n",
    "f = zeros(K); f[50:99] .= 1.0; f[100:139] .= -1.3;\n",
    "X = [x1 zeros(2,K-1)];\n",
    "for k=1:K-1\n",
    "X[:,k+1] = A* X[:,k] + B*f[k]\n",
    "end\n",
    "using Plots\n",
    "plot(X[1,:], xlabel=\"k\", ylabel=\"Position\", legend=false )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/svg+xml": [
       "<?xml version=\"1.0\" encoding=\"utf-8\"?>\n",
       "<svg xmlns=\"http://www.w3.org/2000/svg\" xmlns:xlink=\"http://www.w3.org/1999/xlink\" width=\"600\" height=\"400\" viewBox=\"0 0 2400 1600\">\n",
       "<defs>\n",
       "  <clipPath id=\"clip7000\">\n",
       "    <rect x=\"0\" y=\"0\" width=\"2400\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip7000)\" d=\"\n",
       "M0 1600 L2400 1600 L2400 0 L0 0  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip7001\">\n",
       "    <rect x=\"480\" y=\"0\" width=\"1681\" height=\"1600\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<path clip-path=\"url(#clip7000)\" d=\"\n",
       "M243.864 1425.62 L2352.76 1425.62 L2352.76 47.2441 L243.864 47.2441  Z\n",
       "  \" fill=\"#ffffff\" fill-rule=\"evenodd\" fill-opacity=\"1\"/>\n",
       "<defs>\n",
       "  <clipPath id=\"clip7002\">\n",
       "    <rect x=\"243\" y=\"47\" width=\"2110\" height=\"1379\"/>\n",
       "  </clipPath>\n",
       "</defs>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  300.228,1425.62 300.228,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  632.369,1425.62 632.369,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  964.509,1425.62 964.509,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1296.65,1425.62 1296.65,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1628.79,1425.62 1628.79,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  1960.93,1425.62 1960.93,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  2293.07,1425.62 2293.07,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  243.864,1233.43 2352.76,1233.43 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  243.864,1001.67 2352.76,1001.67 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  243.864,769.916 2352.76,769.916 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  243.864,538.162 2352.76,538.162 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  243.864,306.407 2352.76,306.407 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#000000; stroke-width:2; stroke-opacity:0.1; fill:none\" points=\"\n",
       "  243.864,74.6529 2352.76,74.6529 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  243.864,1425.62 2352.76,1425.62 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  243.864,1425.62 243.864,47.2441 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  300.228,1425.62 300.228,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  632.369,1425.62 632.369,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  964.509,1425.62 964.509,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1296.65,1425.62 1296.65,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1628.79,1425.62 1628.79,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  1960.93,1425.62 1960.93,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  2293.07,1425.62 2293.07,1404.94 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  243.864,1233.43 275.498,1233.43 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  243.864,1001.67 275.498,1001.67 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  243.864,769.916 275.498,769.916 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  243.864,538.162 275.498,538.162 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  243.864,306.407 275.498,306.407 \n",
       "  \"/>\n",
       "<polyline clip-path=\"url(#clip7000)\" style=\"stroke:#000000; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  243.864,74.6529 275.498,74.6529 \n",
       "  \"/>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 300.228, 1479.62)\" x=\"300.228\" y=\"1479.62\">0</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 632.369, 1479.62)\" x=\"632.369\" y=\"1479.62\">100</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 964.509, 1479.62)\" x=\"964.509\" y=\"1479.62\">200</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1296.65, 1479.62)\" x=\"1296.65\" y=\"1479.62\">300</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1628.79, 1479.62)\" x=\"1628.79\" y=\"1479.62\">400</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 1960.93, 1479.62)\" x=\"1960.93\" y=\"1479.62\">500</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:middle;\" transform=\"rotate(0, 2293.07, 1479.62)\" x=\"2293.07\" y=\"1479.62\">600</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 219.864, 1250.93)\" x=\"219.864\" y=\"1250.93\">-0.1</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 219.864, 1019.17)\" x=\"219.864\" y=\"1019.17\">0.0</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 219.864, 787.416)\" x=\"219.864\" y=\"787.416\">0.1</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 219.864, 555.662)\" x=\"219.864\" y=\"555.662\">0.2</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 219.864, 323.907)\" x=\"219.864\" y=\"323.907\">0.3</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:48px; text-anchor:end;\" transform=\"rotate(0, 219.864, 92.1529)\" x=\"219.864\" y=\"92.1529\">0.4</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:66px; text-anchor:middle;\" transform=\"rotate(0, 1298.31, 1559.48)\" x=\"1298.31\" y=\"1559.48\">k</text>\n",
       "</g>\n",
       "<g clip-path=\"url(#clip7000)\">\n",
       "<text style=\"fill:#000000; fill-opacity:1; font-family:Arial,Helvetica Neue,Helvetica,sans-serif; font-size:66px; text-anchor:middle;\" transform=\"rotate(-90, 89.2861, 736.431)\" x=\"89.2861\" y=\"736.431\">Velocity</text>\n",
       "</g>\n",
       "<polyline clip-path=\"url(#clip7002)\" style=\"stroke:#009af9; stroke-width:4; stroke-opacity:1; fill:none\" points=\"\n",
       "  303.55,1001.67 306.871,1001.67 310.193,1001.67 313.514,1001.67 316.835,1001.67 320.157,1001.67 323.478,1001.67 326.8,1001.67 330.121,1001.67 333.442,1001.67 \n",
       "  336.764,1001.67 340.085,1001.67 343.407,1001.67 346.728,1001.67 350.049,1001.67 353.371,1001.67 356.692,1001.67 360.014,1001.67 363.335,1001.67 366.656,1001.67 \n",
       "  369.978,1001.67 373.299,1001.67 376.621,1001.67 379.942,1001.67 383.263,1001.67 386.585,1001.67 389.906,1001.67 393.228,1001.67 396.549,1001.67 399.87,1001.67 \n",
       "  403.192,1001.67 406.513,1001.67 409.835,1001.67 413.156,1001.67 416.477,1001.67 419.799,1001.67 423.12,1001.67 426.442,1001.67 429.763,1001.67 433.085,1001.67 \n",
       "  436.406,1001.67 439.727,1001.67 443.049,1001.67 446.37,1001.67 449.692,1001.67 453.013,1001.67 456.334,1001.67 459.656,1001.67 462.977,1001.67 466.299,1001.67 \n",
       "  469.62,978.495 472.941,955.552 476.263,932.838 479.584,910.35 482.906,888.088 486.227,866.049 489.548,844.229 492.87,822.628 496.191,801.243 499.513,780.072 \n",
       "  502.834,759.113 506.155,738.363 509.477,717.82 512.798,697.483 516.12,677.35 519.441,657.418 522.762,637.685 526.084,618.149 529.405,598.809 532.727,579.662 \n",
       "  536.048,560.707 539.369,541.941 542.691,523.363 546.012,504.97 549.334,486.762 552.655,468.736 555.976,450.889 559.298,433.222 562.619,415.731 565.941,398.415 \n",
       "  569.262,381.272 572.583,364.3 575.905,347.499 579.226,330.865 582.548,314.398 585.869,298.095 589.19,281.955 592.512,265.977 595.833,250.158 599.155,234.498 \n",
       "  602.476,218.994 605.797,203.646 609.119,188.45 612.44,173.407 615.762,158.514 619.083,143.771 622.404,129.174 625.726,114.724 629.047,100.418 632.369,86.2547 \n",
       "  635.69,125.537 639.012,164.426 642.333,202.927 645.654,241.042 648.976,278.777 652.297,316.134 655.619,353.117 658.94,389.731 662.261,425.978 665.583,461.863 \n",
       "  668.904,497.39 672.226,532.56 675.547,567.38 678.868,601.851 682.19,635.977 685.511,669.762 688.833,703.209 692.154,736.322 695.475,769.103 698.797,801.557 \n",
       "  702.118,833.686 705.44,865.494 708.761,896.984 712.082,928.159 715.404,959.022 718.725,989.577 722.047,1019.83 725.368,1049.77 728.689,1079.42 732.011,1108.77 \n",
       "  735.332,1137.83 738.654,1166.59 741.975,1195.07 745.296,1223.27 748.618,1251.18 751.939,1278.81 755.261,1306.17 758.582,1333.25 761.903,1360.06 765.225,1386.61 \n",
       "  768.546,1382.76 771.868,1378.95 775.189,1375.18 778.51,1371.44 781.832,1367.74 785.153,1364.08 788.475,1360.46 791.796,1356.87 795.117,1353.32 798.439,1349.8 \n",
       "  801.76,1346.32 805.082,1342.87 808.403,1339.46 811.724,1336.08 815.046,1332.74 818.367,1329.43 821.689,1326.15 825.01,1322.91 828.331,1319.69 831.653,1316.51 \n",
       "  834.974,1313.37 838.296,1310.25 841.617,1307.16 844.939,1304.11 848.26,1301.08 851.581,1298.09 854.903,1295.12 858.224,1292.19 861.546,1289.29 864.867,1286.41 \n",
       "  868.188,1283.56 871.51,1280.74 874.831,1277.95 878.153,1275.19 881.474,1272.45 884.795,1269.75 888.117,1267.07 891.438,1264.41 894.76,1261.78 898.081,1259.18 \n",
       "  901.402,1256.61 904.724,1254.06 908.045,1251.53 911.367,1249.04 914.688,1246.56 918.009,1244.11 921.331,1241.69 924.652,1239.29 927.974,1236.91 931.295,1234.56 \n",
       "  934.616,1232.23 937.938,1229.93 941.259,1227.64 944.581,1225.38 947.902,1223.15 951.223,1220.93 954.545,1218.74 957.866,1216.57 961.188,1214.42 964.509,1212.29 \n",
       "  967.83,1210.19 971.152,1208.1 974.473,1206.04 977.795,1203.99 981.116,1201.97 984.437,1199.97 987.759,1197.98 991.08,1196.02 994.402,1194.08 997.723,1192.15 \n",
       "  1001.04,1190.25 1004.37,1188.36 1007.69,1186.5 1011.01,1184.65 1014.33,1182.82 1017.65,1181.01 1020.97,1179.21 1024.29,1177.44 1027.62,1175.68 1030.94,1173.94 \n",
       "  1034.26,1172.22 1037.58,1170.51 1040.9,1168.82 1044.22,1167.15 1047.54,1165.5 1050.87,1163.86 1054.19,1162.24 1057.51,1160.63 1060.83,1159.04 1064.15,1157.47 \n",
       "  1067.47,1155.91 1070.79,1154.37 1074.12,1152.84 1077.44,1151.33 1080.76,1149.83 1084.08,1148.35 1087.4,1146.88 1090.72,1145.43 1094.04,1143.99 1097.37,1142.57 \n",
       "  1100.69,1141.16 1104.01,1139.77 1107.33,1138.39 1110.65,1137.02 1113.97,1135.66 1117.29,1134.32 1120.61,1133 1123.94,1131.69 1127.26,1130.38 1130.58,1129.1 \n",
       "  1133.9,1127.82 1137.22,1126.56 1140.54,1125.31 1143.86,1124.08 1147.19,1122.85 1150.51,1121.64 1153.83,1120.44 1157.15,1119.25 1160.47,1118.08 1163.79,1116.91 \n",
       "  1167.11,1115.76 1170.44,1114.62 1173.76,1113.49 1177.08,1112.37 1180.4,1111.27 1183.72,1110.17 1187.04,1109.08 1190.36,1108.01 1193.69,1106.95 1197.01,1105.89 \n",
       "  1200.33,1104.85 1203.65,1103.82 1206.97,1102.8 1210.29,1101.79 1213.61,1100.79 1216.94,1099.8 1220.26,1098.81 1223.58,1097.84 1226.9,1096.88 1230.22,1095.93 \n",
       "  1233.54,1094.99 1236.86,1094.05 1240.19,1093.13 1243.51,1092.21 1246.83,1091.31 1250.15,1090.41 1253.47,1089.53 1256.79,1088.65 1260.11,1087.78 1263.44,1086.92 \n",
       "  1266.76,1086.06 1270.08,1085.22 1273.4,1084.38 1276.72,1083.56 1280.04,1082.74 1283.36,1081.93 1286.69,1081.12 1290.01,1080.33 1293.33,1079.54 1296.65,1078.76 \n",
       "  1299.97,1077.99 1303.29,1077.23 1306.61,1076.48 1309.93,1075.73 1313.26,1074.99 1316.58,1074.25 1319.9,1073.53 1323.22,1072.81 1326.54,1072.1 1329.86,1071.39 \n",
       "  1333.18,1070.7 1336.51,1070.01 1339.83,1069.32 1343.15,1068.65 1346.47,1067.98 1349.79,1067.31 1353.11,1066.66 1356.43,1066.01 1359.76,1065.36 1363.08,1064.73 \n",
       "  1366.4,1064.1 1369.72,1063.47 1373.04,1062.85 1376.36,1062.24 1379.68,1061.64 1383.01,1061.04 1386.33,1060.44 1389.65,1059.86 1392.97,1059.27 1396.29,1058.7 \n",
       "  1399.61,1058.13 1402.93,1057.56 1406.26,1057 1409.58,1056.45 1412.9,1055.9 1416.22,1055.36 1419.54,1054.82 1422.86,1054.29 1426.18,1053.77 1429.51,1053.24 \n",
       "  1432.83,1052.73 1436.15,1052.22 1439.47,1051.71 1442.79,1051.21 1446.11,1050.72 1449.43,1050.23 1452.76,1049.74 1456.08,1049.26 1459.4,1048.78 1462.72,1048.31 \n",
       "  1466.04,1047.85 1469.36,1047.39 1472.68,1046.93 1476.01,1046.48 1479.33,1046.03 1482.65,1045.58 1485.97,1045.14 1489.29,1044.71 1492.61,1044.28 1495.93,1043.85 \n",
       "  1499.25,1043.43 1502.58,1043.01 1505.9,1042.6 1509.22,1042.19 1512.54,1041.79 1515.86,1041.38 1519.18,1040.99 1522.5,1040.59 1525.83,1040.21 1529.15,1039.82 \n",
       "  1532.47,1039.44 1535.79,1039.06 1539.11,1038.69 1542.43,1038.32 1545.75,1037.95 1549.08,1037.59 1552.4,1037.23 1555.72,1036.87 1559.04,1036.52 1562.36,1036.17 \n",
       "  1565.68,1035.83 1569,1035.49 1572.33,1035.15 1575.65,1034.81 1578.97,1034.48 1582.29,1034.15 1585.61,1033.83 1588.93,1033.51 1592.25,1033.19 1595.58,1032.87 \n",
       "  1598.9,1032.56 1602.22,1032.25 1605.54,1031.95 1608.86,1031.64 1612.18,1031.34 1615.5,1031.05 1618.83,1030.75 1622.15,1030.46 1625.47,1030.17 1628.79,1029.89 \n",
       "  1632.11,1029.61 1635.43,1029.33 1638.75,1029.05 1642.08,1028.78 1645.4,1028.51 1648.72,1028.24 1652.04,1027.97 1655.36,1027.71 1658.68,1027.45 1662,1027.19 \n",
       "  1665.33,1026.94 1668.65,1026.68 1671.97,1026.43 1675.29,1026.19 1678.61,1025.94 1681.93,1025.7 1685.25,1025.46 1688.57,1025.22 1691.9,1024.98 1695.22,1024.75 \n",
       "  1698.54,1024.52 1701.86,1024.29 1705.18,1024.07 1708.5,1023.84 1711.82,1023.62 1715.15,1023.4 1718.47,1023.18 1721.79,1022.97 1725.11,1022.76 1728.43,1022.54 \n",
       "  1731.75,1022.34 1735.07,1022.13 1738.4,1021.92 1741.72,1021.72 1745.04,1021.52 1748.36,1021.32 1751.68,1021.13 1755,1020.93 1758.32,1020.74 1761.65,1020.55 \n",
       "  1764.97,1020.36 1768.29,1020.17 1771.61,1019.99 1774.93,1019.8 1778.25,1019.62 1781.57,1019.44 1784.9,1019.27 1788.22,1019.09 1791.54,1018.92 1794.86,1018.74 \n",
       "  1798.18,1018.57 1801.5,1018.4 1804.82,1018.24 1808.15,1018.07 1811.47,1017.91 1814.79,1017.74 1818.11,1017.58 1821.43,1017.42 1824.75,1017.27 1828.07,1017.11 \n",
       "  1831.4,1016.96 1834.72,1016.8 1838.04,1016.65 1841.36,1016.5 1844.68,1016.35 1848,1016.21 1851.32,1016.06 1854.65,1015.92 1857.97,1015.78 1861.29,1015.63 \n",
       "  1864.61,1015.5 1867.93,1015.36 1871.25,1015.22 1874.57,1015.08 1877.89,1014.95 1881.22,1014.82 1884.54,1014.69 1887.86,1014.56 1891.18,1014.43 1894.5,1014.3 \n",
       "  1897.82,1014.17 1901.14,1014.05 1904.47,1013.92 1907.79,1013.8 1911.11,1013.68 1914.43,1013.56 1917.75,1013.44 1921.07,1013.32 1924.39,1013.21 1927.72,1013.09 \n",
       "  1931.04,1012.98 1934.36,1012.86 1937.68,1012.75 1941,1012.64 1944.32,1012.53 1947.64,1012.42 1950.97,1012.32 1954.29,1012.21 1957.61,1012.1 1960.93,1012 \n",
       "  1964.25,1011.9 1967.57,1011.79 1970.89,1011.69 1974.22,1011.59 1977.54,1011.49 1980.86,1011.4 1984.18,1011.3 1987.5,1011.2 1990.82,1011.11 1994.14,1011.01 \n",
       "  1997.47,1010.92 2000.79,1010.83 2004.11,1010.73 2007.43,1010.64 2010.75,1010.55 2014.07,1010.47 2017.39,1010.38 2020.72,1010.29 2024.04,1010.2 2027.36,1010.12 \n",
       "  2030.68,1010.03 2034,1009.95 2037.32,1009.87 2040.64,1009.79 2043.97,1009.71 2047.29,1009.62 2050.61,1009.55 2053.93,1009.47 2057.25,1009.39 2060.57,1009.31 \n",
       "  2063.89,1009.23 2067.21,1009.16 2070.54,1009.08 2073.86,1009.01 2077.18,1008.94 2080.5,1008.86 2083.82,1008.79 2087.14,1008.72 2090.46,1008.65 2093.79,1008.58 \n",
       "  2097.11,1008.51 2100.43,1008.44 2103.75,1008.38 2107.07,1008.31 2110.39,1008.24 2113.71,1008.18 2117.04,1008.11 2120.36,1008.05 2123.68,1007.98 2127,1007.92 \n",
       "  2130.32,1007.86 2133.64,1007.8 2136.96,1007.73 2140.29,1007.67 2143.61,1007.61 2146.93,1007.55 2150.25,1007.5 2153.57,1007.44 2156.89,1007.38 2160.21,1007.32 \n",
       "  2163.54,1007.27 2166.86,1007.21 2170.18,1007.15 2173.5,1007.1 2176.82,1007.05 2180.14,1006.99 2183.46,1006.94 2186.79,1006.89 2190.11,1006.83 2193.43,1006.78 \n",
       "  2196.75,1006.73 2200.07,1006.68 2203.39,1006.63 2206.71,1006.58 2210.04,1006.53 2213.36,1006.48 2216.68,1006.43 2220,1006.39 2223.32,1006.34 2226.64,1006.29 \n",
       "  2229.96,1006.25 2233.29,1006.2 2236.61,1006.16 2239.93,1006.11 2243.25,1006.07 2246.57,1006.02 2249.89,1005.98 2253.21,1005.94 2256.53,1005.89 2259.86,1005.85 \n",
       "  2263.18,1005.81 2266.5,1005.77 2269.82,1005.73 2273.14,1005.69 2276.46,1005.65 2279.78,1005.61 2283.11,1005.57 2286.43,1005.53 2289.75,1005.49 2293.07,1005.45 \n",
       "  \n",
       "  \"/>\n",
       "</svg>\n"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "plot(X[2,:], xlabel=\"k\", ylabel=\"Velocity\", legend=false )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Figure 9.4** Simulation of a mass moving along a line: position (top) and\n",
    "velocity (bottom)."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.1.1",
   "language": "julia",
   "name": "julia-1.1"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.1.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
